Introduction à la modélisation statistique bayésienne

Ce cours est grandement inspiré des livres suivants :

• McElreath, R. (2016, 2020). Statistical Rethinking: A Bayesian Course with Examples in R and Stan. CRC Press.

• Kurz, S. (2019). Statistical Rethinking with brms, ggplot2, and the tidyverse. Available online.

• Kruschke, J. K. (2015). Doing Bayesian Data Analysis, Second Edition: A Tutorial with R, JAGS, and Stan. Academic Press / Elsevier.

• Gelman, A., Carlin, J. B., Stern, H. S., Dunson, D. B., Vehtari, A., & Rubin, D. B. (2013). Bayesian Data Analysis, third edition. London: CRC Press.

• Lambert, B. (2018). A Student’s Guide to Bayesian Statistics. SAGE Publications Ltd.

• Nicenboim, B., Schad, D., & Vasishth, S. (2021). An Introduction to Bayesian Data Analysis for Cognitive Science. Available online.

Les slides, données, et scripts seront disponibles juste avant chaque séance sur Github : https://github.com/lnalborczyk/IMSB2021

Objectifs généraux

• Comprendre les concepts fondamentaux de la statistique bayésienne
• Être capable d'implémenter des modèles bayésiens simples en R
• Réaliser que l'approche bayésienne est plus intuitive que l'approche fréquentiste

Objectifs pratiques

• Être capable de réaliser une analyse complète (i.e., identification du modèle approprié, écriture du modèle mathématique, implémentation en R, interprétation et report des résultats) d'un jeu de données simple

Cours n°01 : Introduction à l'inférence bayésienne
Cours n°02 : Modèle Beta-Binomial
Cours n°03 : Introduction à brms, modèle de régression linéaire
Cours n°04 : Modèle de régression linéaire (suite)
Cours n°05 : Markov Chain Monte Carlo
Cours n°06 : Modèle linéaire généralisé
Cours n°07 : Comparaison de modèles
Cours n°08 : Modèles multi-niveaux
Cours n°09 : Modèles multi-niveaux généralisés
Cours n°10 : Data Hackathon

Une probabilité est une valeur numérique assignée à un événement \( A \), compris comme une possibilité appartenant à l'univers \( \Omega \) (l'ensemble de toutes les issues possibles). Les probabilités se conforment aux axiomes suivants :

• Non-negativity : \( \Pr(A) \geq 0 \)
• Normalisation : \( \Pr(\Omega) = 1 \)
• Additivity (pour des événements incompatibles) : \( \Pr(A_{1} \cup A_{2}) = \Pr(A_{1}) + \Pr(A_{2}) \)

Le dernier axiome est également connu comme la règle de la somme, et peut se généraliser à des événements non mutuellement exclusifs : \( \Pr(A_{1} \cup A_{2}) = \Pr(A_{1}) + \Pr(A_{2}) - \Pr(A_{1} \cap A_{2}) \).

Quelle est la probabilité…

• D'obtenir un chiffre pair sur un lancer de dé ?
• Que j'apprenne quelque chose pendant cette formation ?

Est-ce qu'il s'agit, pour chaque exemple, de la même sorte de probabilité ?

\[ \Pr(\text{pair}) = \frac{\text{nombre de cas favorables}}{\text{nombre de cas possibles}} = \frac{3}{6} = \frac{1}{2} \]

Problème : cette définition est uniquement applicable aux situations dans lesquelles il n'y a qu'un nombre fini de résultats possibles équiprobables…

Quelle est la probabilité qu'il pleuve demain ?

\[ \Pr(\text{pluie}) = \frac{\text{pluie}}{ \{ \text{pluie, non-pluie} \} } = \frac{1}{2} \]

\[ \Pr(x) = \lim_{n_{t} \to \infty}\frac{n_{x}}{n_{t}} \]

Où \( n_{x} \) est le nombre d'occurrences de l’événement \( x \) et \( n_{t} \) le nombre total d'essais. L'interprétation fréquentiste postule que, à long-terme (i.e., quand le nombre d'essais s'approche de l'infini), la fréquence relative va converger exactement vers ce qu'on appelle “probabilité”.

Conséquence : le concept de probabilité s'applique uniquement aux collectifs, et non aux événements singuliers.

library(tidyverse)

sample(c(0, 1), 500, replace = TRUE) %>%
        data.frame %>%
        mutate(x = seq_along(.), y = cumsum(.) / seq_along(.) ) %>%
        ggplot(aes(x = x, y = y), log = "y") +
        geom_line(lwd = 1) +
        geom_hline(yintercept = 0.5, lty = 3) +
        xlab("Nombre de lancers") +
        ylab("Proportion de faces") +
        ylim(0, 1)

Quelle classe de référence ? Quelle est la probabilité que je vive jusqu'à 80 ans ? En tant qu'homme ? En tant que Français ?

Quid des événements qui ne peuvent pas se répéter ? Quelle est la probabilité que j'apprenne quelque chose pendant cette formation ?

À partir de combien de lancers (d'une pièce par exemple) a-t-on une bonne approximation de la probabilité ? Une classe finie d’événements de taille \( n \) ne peut produire que des fréquences relatives de précision \( 1/n \)…

Les propriétés fréquentistes (i.e., à long terme) des objets (e.g., une pièce) seraient provoquées par des propriétés physiques intrinsèques aux objets. Par exemple, une pièce biaisée va engendrer une fréquence relative (et donc une probabilité) biaisée en raison de ses propriétés physiques. Pour les propensionnistes, les probabilités représentent ces caractéristiques intrinsèques, ces propensions à générer certaines fréquences relatives, et non les fréquences relatives en elles-mêmes.

Conséquence : ces propriétés sont les propriétés d’événements individuels… et non de séquences ! L’interprétation propensionniste nous permet donc de parler de la probabilité d'événements uniques.

Conclusion n°1 : l'étudiant tiré au sort porte un t-shirt ✔

Conclusion n°2 : l'étudiant tiré au sort porte un t-shirt vert ✗

Conclusion n°3 : l'étudiant tiré au sort porte un t-shirt rouge ✗

L'interprétation logique du concept de probabilité essaye de généraliser la logique (vrai / faux) au monde probabiliste. La probabilité représente donc le degré de support logique qu'une conclusion peut avoir, relativement à un ensemble de prémisses (Keynes, 1921; Carnap, 1950).

Conséquence : toute probabilité est conditionnelle.

La probabilité est une mesure du degré d'incertitude. Un événement certain aura donc une probabilité de 1 et un événement impossible aura une probabilité de 0.

Pour parler de probabilités, dans ce cadre, nous n'avons donc plus besoin de nous référer à la limite d'occurrence d'un événement (fréquence).

So to assign equal probabilities to two events is not in any way an assertion that they must occur equally often in any random experiment […], it is only a formal way of saying I don’t know (Jaynes, 1986).

• Interprétation classique (Laplace, Bernouilli, Leibniz)
• Interprétation fréquentiste (Venn, Reichenbach, von Mises)
• Interprétation propensionniste (Popper, Miller)
• Interprétation logique (Keynes, Carnap)
• Interprétation bayésienne (Jeffreys, de Finetti, Savage)

Plus de détails sur la Stanford Encyclopedia of Philosophy…

Probabilité épistémique

• Toute probabilité est conditionnelle à de l'information disponible (e.g., prémisses ou données)
• La probabilité comme moyen de quantifier l'incertitude
• Probabilité logique, probabilité bayésienne

Probabilité physique

• Les probabilités dépendent d'un état du monde, de caractéristiques physiques, elles sont indépendantes de l'information disponible (ou de l'incertitude)
• Probabilité classique, probabilité fréquentiste

Une suite de nombres au hasard…

[1] 79 27 98 28 42

runif, rnorm, rbinom...

RNGkind()[1] # default random number generator

[1] "Mersenne-Twister"

Inventé par Makoto Matsumoto et Takuji Nishimura en 1997, le Mersenne Twister est un générateur de nombres pseudo-aléatoires. L'algorithme Mersenne-Twister est (était) le générateur par défaut en Python, Ruby, R, Matlab…

\[ Y_{i} \sim \mathrm{Uniform}(10,100) \]

set.seed(666)
as.integer(runif(5, 10, 100) )

[1] 79 27 98 28 42

Si je connais la graine d'origine (i.e., le “départ” de l'algorithme) et le fonctionnement précis de l'algorithme, je peux prédire les nombres qui vont être générés. Il sera alors difficile de soutenir que ces nombres auront été générés “au hasard”…

[1] 79 27 98 28 42 76

[1] 79 27 98 28 42 76 98

“Le hasard est un concept indiquant l'aspect aléatoire des choses, des événements ou encore des décisions, signalant ainsi l'impossibilité de prévoir avec certitude un fait quelconque. Ainsi, pour éclairer le sens du mot, il peut-être dit que le hasard est synonyme d'« imprévisibilité », ou « imprédictibilité »” (Wikipédia).

Ce que nous appelons hasard n’est et ne peut être que la cause ignorée d’un effet connu. – Voltaire

Le hasard, ce sont les lois que nous ne connaissons pas. – Émile Borel

Randomness is a proxy for lack of knowledge. – Richard McElreath

Le hasard apparaît donc comme un état subjectif d'indétermination des causes. Pour parler de l'incertain (i.e., pour quantifier l'incertitude), on utilise les probabilités comme langage.

Exemple 1

• Si un suspect ment, il transpire. (On observe que) Ce suspect transpire.
• Par conséquent, ce suspect ment.

Exemple 2

• Si un suspect transpire, il ment. (On observe que) Ce suspect ne transpire pas.
• Par conséquent, ce suspect ne ment pas.

Exemple 3

• Tous les menteurs transpirent. (On observe que) Ce suspect ne transpire pas.
• Par conséquent, ce suspect n'est pas un menteur.

• Affirmation du conséquent : \( \dfrac{A \Rightarrow B, \ B}{A} \)

• Si il a plu, alors le sol est mouillé (A implique B). Le sol est mouillé (B). Donc il a plu (A).

• Négation de l'antécédent : \( \dfrac{A \Rightarrow B, \ \neg A}{\neg B} \)

• Si il a plu, alors le sol est mouillé (A implique B). Il n'a pas plus (non A). Donc le sol n'est pas mouillé (non B).

• Modus ponens : \( \dfrac{A \Rightarrow B, \ A}{B} \)

• Si on est lundi, alors John ira au travail (A implique B). On est lundi (A). Donc John ira au travail (B).

• Modus tollens : \( \dfrac{A \Rightarrow B, \ \neg B}{\neg A} \)

• Si mon chien détecte un intru, alors il aboie (A implique B). Mon chien n'a pas aboyé (non B). Donc il n'a pas détecté d'intrus (non A).

Le modus tollens est un des raisonnements logiques les plus importants et les plus performants. Dans le cadre de l'inférence statistique, il s'applique parfaitement au cas suivant : “Si \( H_{0} \) est vraie, alors \( x \) ne devrait pas se produire. On observe \( x \). Alors \( H_{0} \) est fausse”.

Mais nous avons le plus souvent affaire à des hypothèses “continues”, probabilistes.

L'inférence fréquentiste (Fishérienne) est elle aussi probabiliste, de la forme “Si \( H_{0} \) est vraie, alors \( x \) est peu probable. On observe \( x \). Alors \( H_{0} \) est peu probable.”

Or cet argument est invalide, le modus tollens ne s'applique pas au monde probabiliste (e.g., Pollard & Richardson, 1987; Rouder, Morey, Verhagen, Province, & Wagenmakers, 2016).

Par exemple : si un individu est un homme, alors il est peu probable qu'il soit pape. François est pape. François n'est donc certainement pas un homme…

Poppérisme naïf : la science progresse par falsification logique, donc la statistique devrait viser la falsification. Mais…

• Les hypothèses théoriques ne sont pas les modèles (hypothèses statistiques)

Models are devices that connect theories to data. A model is an instanciation of a theory as a set of probabilistic statements (Rouder et al., 2016).

• Nos hypothèses sont souvent probabilistes
• Erreurs de mesure
• La falsification concerne le problème de la démarcation, pas celui de la méthode
• La science est une technologie sociale, la falsification est consensuelle, et non pas logique

On s'intéresse au lien entre deux variables aléatoires continues, \( x \) et \( y \).

L'hypothèse de modélisation la plus simple est de postuler une relation linéaire.

Cette description peut-être améliorée pour mieux prendre en compte les données qui s'écartent de la prédiction linéaire.

Un ensemble de \( N \) points peut être exhaustivement (i.e., sans erreur) décrit par une fonction polynomiale d'ordre \( N-1 \). Augmenter la complexité du modèle améliore donc la précision de notre description des données mais réduit la généralisabilité de ses prédictions (bias-variance tradeoff).

Nous avons besoin d'outils qui prennent en compte le rapport qualité de la description / complexité, c'est à dire la parcimonie des modèles (e.g., AIC, WAIC).

Besoin d'un cadre pour développer des modèles cohérents. Nos outils :

• Bayesian data analysis, utiliser les probabilités pour décrire l'incertitude
• Multilevel modelling, des modèles à multiples niveaux d'incertitude
• Approche par comparaison de modèles, au lieu d'essayer de falsifier un null model, comparer des modèles intéressants (AIC, WAIC)

Imaginons que nous disposions d'un sac contenant 4 billes. Ces billes peuvent être soit blanches, soit bleues. Nous savons qu'il y a précisément 4 billes, mais nous ne connaissons pas le nombre de billes de chaque couleur.

Nous savons cependant qu'il existe cinq possibilités (que nous considérons comme nos hypothèses) :

⚪ ⚪ ⚪ ⚪

🔵 ⚪ ⚪ ⚪

🔵 🔵 ⚪ ⚪

🔵 🔵 🔵 ⚪

🔵 🔵 🔵 🔵

Le but est de déterminer quelle combinaison serait la plus probable, sachant certaines observations. Imaginons que l'on tire trois billes à la suite, avec remise, et que l'on obtienne la séquence suivante :

🔵 ⚪ 🔵

Cette séquence représente nos données. À partir de là, quelle inférence peut-on faire sur le contenu du sac ? En d'autres termes, que peut-on dire de la probabilité de chaque hypothèse ?

⚪ ⚪ ⚪ ⚪

🔵 ⚪ ⚪ ⚪

🔵 🔵 ⚪ ⚪

🔵 🔵 🔵 ⚪

🔵 🔵 🔵 🔵

hypothèse: 🔵 ⚪ ⚪ ⚪                  données: 🔵

hypothèse : 🔵 ⚪ ⚪ ⚪                  données : 🔵 ⚪

hypothèse : 🔵 ⚪ ⚪ ⚪                  données : 🔵 ⚪ 🔵

hypothèse : 🔵 ⚪ ⚪ ⚪                  données : 🔵 ⚪ 🔵

Sous cette hypothèse, \( 3 \) chemins (sur \( 4^{3}=64 \)) conduisent au résultat obtenu. Qu'en est-il des autres hypothèses?

⚪ ⚪ ⚪ 🔵      ⚪ 🔵 🔵 🔵      ⚪ ⚪ 🔵 🔵

🔵 🔵 🔵 ⚪

Au vu des données, cette hypothèse est la plus probable car c'est l'hypothèse qui maximise le nombre de manières possibles d'obtenir les données obtenues.

Dans le cas précédent, nous avons considéré que toutes les hypothèses étaient équiprobables a priori (suivant le principe d'indifférence). Cependant, on pourrait avoir de l'information a priori, provenant de nos connaissances (des particularités des sacs de billes par exemple) ou de données antérieures.

Imaginons que nous tirions une nouvelle bille du sac, comment incorporer cette nouvelle donnée ?

Il suffit d'appliquer la même stratégie que précédemment, et de mettre à jour le dernier compte en le multipliant par ces nouvelles données. Yesterday's posterior is today's prior (Lindley, 2000).

Supposons maintenant qu'un employé de l'usine de fabrication des billes nous dise que les billes bleues sont rares… Cet employé nous dit que pour chaque sac contenant 3 billes bleues, ils fabriquent deux sacs en contenant seulement deux, et trois sacs en contenant seulement une. Il nous apprend également que tous les sacs contiennent au moins une bille bleue et une bille blanche…

La probabilité d'une hypothèse après avoir observé certaines données est proportionnelle au nombre de façons qu'a cette hypothèse de produire les données observées, multiplié par sa probabilité a priori.

\[ \Pr(\text{hypothèse}\ |\ \text{données}) \propto \Pr(\text{données}\ |\ \text{hypothèse}) \times \Pr(\text{hypothèse}) \]

Pour passer des plausibilités aux probabilités, il suffit de standardiser ces plausibilités pour que la somme des plausibilités de toutes les hypothèses possibles soit égale à \( 1 \).

\[ \Pr(\text{hypothèse}\ |\ \text{données}) = \frac{\Pr(\text{données}\ |\ \text{hypothèse})\times \Pr(\text{hypothèse})}{\text{somme des produits}} \]

Définissons \( p \) comme la proportion de billes bleues.

ways <- c(0, 3, 8, 9, 0)
ways / sum(ways)

[1] 0.00 0.15 0.40 0.45 0.00

• \( \theta \) est un paramètre ou vecteur de paramètres (e.g., la proportion de billes bleues).
• \( \color{orangered}{p(x\vert \theta)} \) la probabilité conditionnelle des données \( x \) sachant le paramètre \( \theta \) \( \color{orangered}{[p(x | \theta = \theta)]} \).
• \( \color{orangered}{p(x\vert \theta)} \) une fois que la valeur de \( x \) est connue, est vue comme la fonction de vraisemblance (likelihood) du paramètre \( \theta \). Attention, il ne s'agit pas d'une distribution de probabilité \( \color{orangered}{[p(x = x | \theta)]} \).
• \( \color{steelblue}{p(\theta)} \) la probabilité a priori de \( \theta \).
• \( \color{purple}{p(\theta \vert x)} \) la probabilité a posteriori de \( \theta \) (sachant \( x \)).
• \( \color{green}{p(x)} \) la probabilité marginale de \( x \) (sur \( \theta \)) ou “vraisemblance marginale”, “vraisemblance intégrée”.

\[ \color{purple}{p(\theta \vert x)} = \dfrac{\color{orangered}{p(x\vert \theta)} \color{steelblue}{p(\theta)}}{\color{green}{p(x)}} = \dfrac{\color{orangered}{p(x\vert \theta)} \color{steelblue}{p(\theta)}}{\color{green}{\sum\limits_{\theta}p(x|\theta)p(\theta)}} = \dfrac{\color{orangered}{p(x\vert \theta)} \color{steelblue}{p(\theta)}}{\color{green}{\int\limits_{\theta}p(x|\theta)p(\theta)\mathrm{d}x}} \propto \color{orangered}{p(x\vert \theta)} \color{steelblue}{p(\theta)} \]

Dans ce cadre, pour chaque problème, nous allons suivre 3 étapes :

• Construire le modèle (l'histoire des données): likelihood + prior
• Mettre à jour grâce aux données (updating), calculer la probabilité a posteriori
• Evaluer le modèle, fit, sensibilité, résumer les résultats, ré-ajuster

Bayesian inference is really just counting and comparing of possibilities […] in order to make good inference about what actually happened, it helps to consider everything that could have happened (McElreath, 2015).

Une fonction de masse (probability mass function, ou PMF) est une fonction qui attribue une probabilité à chaque valeur d'une variable aléatoire. Exemple de la distribution binomiale pour une pièce non biaisée (\( \theta = 0.5 \)), probabilité d'obtenir \( N \) faces sur 10 lancers.

# PMFs sum to 1
dbinom(x = 0:10, size = 10, prob = 0.5) %>% sum

[1] 1

Une (fonction de) densité de probabilité (probability density function, ou PDF), est une fonction qui permet de représenter une loi de probabilité sous forme d'intégrales (l'équivalent de la PMF pour des variables aléatoires strictement continues).

# PDFs integrate to 1
integrate(dnorm, -Inf, Inf, mean = 100, sd = 15)

1 with absolute error < 1.3e-06

Une intégrale correspond à la surface (aire géométrique) délimitée par la représentation graphique d'une fonction, l'aire sous la courbe. Une distribution est dite impropre si son intégrale n'est pas égale à un nombre fini (e.g., \( + \infty \)) et normalisée si son intégrale est égale à 1.

L'intégrale de \( f(x) \) sur l'intervalle [90 ; 96] vaut : \( p(90 < x < 96) = \int_{90}^{96} f(x) \ \mathrm{d}x = 0.142 \).

integrate(dnorm, 90, 96, mean = 100, sd = 15)

0.1423704 with absolute error < 1.6e-15

library(tidyverse)

data(HairEyeColor) # data adapted from Snee (1974)

cont <- apply(HairEyeColor, c(1, 2), sum) %>% t 
cont <- round(cont / sum(cont), 2)
cont

       Hair
Eye     Black Brown  Red Blond
  Brown  0.11  0.20 0.04  0.01
  Blue   0.03  0.14 0.03  0.16
  Hazel  0.03  0.09 0.02  0.02
  Green  0.01  0.05 0.02  0.03

Dans chaque cellule, on trouve la probabilité conjointe d'avoir telle couleur de cheveux ET telle couleur d'yeux, qui s'écrit \( p(c,y)=p(y,c) \).

cont2 <- cont %>% as.data.frame %>% mutate(marginal_eye = rowSums(cont) )
rownames(cont2) <- row.names(cont)
cont2

      Black Brown  Red Blond marginal_eye
Brown  0.11  0.20 0.04  0.01         0.36
Blue   0.03  0.14 0.03  0.16         0.36
Hazel  0.03  0.09 0.02  0.02         0.16
Green  0.01  0.05 0.02  0.03         0.11

On peut aussi s'intéresser à la probabilité d'avoir des yeux bleus, de manière générale. Il s'agit de la probabilité marginale de l'événement yeux bleus, qui s'obtient par la somme de toutes les probabilités jointes impliquant l'événement yeux bleus. Elle s'écrit \( p(y)=\sum\limits_{c}p(y|c)p(c) \).

cont3 <- rbind(cont2, colSums(cont2) )
rownames(cont3) <- c(row.names(cont2), "marginal_hair")
cont3

              Black Brown  Red Blond marginal_eye
Brown          0.11  0.20 0.04  0.01         0.36
Blue           0.03  0.14 0.03  0.16         0.36
Hazel          0.03  0.09 0.02  0.02         0.16
Green          0.01  0.05 0.02  0.03         0.11
marginal_hair  0.18  0.48 0.11  0.22         0.99

On peut bien entendu aussi s'intéresser aux probabilités des couleurs de cheveux, de manière générale. Elle s'écrit \( p(c)=\sum\limits_{y}p(c|y)p(y) \).

On pourrait aussi s'intéresser à la probabilité qu'une personne ait les cheveux blonds, sachant qu'elle a les yeux bleus. Il s'agit d'une probabilité conditionnelle, et s'écrit \( p(c|y) \). Cette probabilité conditionnelle peut se ré-écrire : \( p(c | y) = \frac{p(c, y)}{p(y)} \).

     Black Brown  Red Blond marginal_eye
Blue  0.03  0.14 0.03  0.16         0.36

Par exemple, quelle est la probabilité d'avoir des yeux bleus lorsqu'on a les cheveux blonds ?

cont3["Blue", "Blond"] / cont3["Blue", "marginal_eye"]  

     Blue 
0.4444444 

On remarque dans le cas précédent que \( p(blonds|bleus) \) n'est pas nécessairement égal à \( p(bleus|blonds) \).

Autre exemple : la probabilité de mourir sachant qu'on a été attaqué par un requin n'est pas la même que la probabilité d'avoir été attaqué par un requin, sachant qu'on est mort (confusion of the inverse). De la même manière, \( p(data|H_{0}) \neq p(H_{0}|data) \).

À partir des axiomes de Kolmogorov (cf. début du cours), et des définitions précédentes des probabilités conjointes, marginales, et conditionnelles, découle la règle du produit :

\[ p(a, b) = p(b) \cdot p(a|b) = p(a) \cdot p(b|a) \]

\[ p(x,y) = p(x|y)p(y) = p(y|x)p(x) \]

\[ p(y|x)p(x) = p(x|y)p(y) \]

\[ p(y|x) = \dfrac{p(x|y)p(y)}{p(x)} \]

\[ p(x|y) = \dfrac{p(y|x)p(x)}{p(y)} \]

\[ \Pr(\text{hypothèse}\ |\ \text{données}) = \frac{\Pr(\text{données}\ |\ \text{hypothèse})\times \Pr(\text{hypothèse})}{\text{somme des produits}} \]

• Chez les femmes âgées de 40-50 ans, sans antécédents familiaux et sans symptômes, la probabilité d'avoir un cancer du sein est de 0.008.

• Propriétés de la mammographie :

  • Si une femme a un cancer du sein, la probabilité d'avoir un résultat positif est de 0.90.
  • Si une femme n'a pas de cancer du sein, la probabilité d'avoir un résultat positif est de 0.07.

• Imaginons qu'une femme passe une mammographie, et que le test est positif. Que doit-on inférer ? Quelle est la probabilité que cette femme ait un cancer du sein ?

• Une approche générale de l'estimation de paramètre
• Les paramètres gouvernent les données, les données dépendent des paramètres
  • Sachant certaines valeurs des paramètres, nous pouvons calculer la probabilité conditionnelle des données observées
  • Le résultat de la mammographie (i.e., les données) dépend de la présence / absence d'un cancer du sein (i.e., le paramètre)
• L'approche par maximum de vraisemblance pose la question : “Quelles sont les valeurs du paramètre qui rendent les données observées les plus probables ?”
• Spécifier la probabilité conditionnelle des données \( p(x|\theta) \)
• Quand on le considère comme fonction de \( \theta \), on parle de likelihood : \( L(\theta|x) = p(X = x|\theta) \)
• L'approche par maximum de vraisemblance consiste donc à maximiser cette fonction, en utilisant les valeurs (connues) de \( x \)

• Si une femme a un cancer du sein, la probabilité d'obtenir un résultat positif est de .90

  • \( \Pr(\text{Mam=+} | \text{Cancer=+}) = 0.90 \)
  • \( \Pr(\text{Mam=-} | \text{Cancer=+}) = 0.10 \)

• Si une femme n'a pas de cancer du sein, la probabilité d'obtenir un résultat positif est de .07

  • \( \Pr(\text{Mam=+} | \text{Cancer=-}) = 0.07 \)
  • \( \Pr(\text{Mam=-} | \text{Cancer=-}) = 0.93 \)

• Une femme passe une mammographie, le résultat est positif…

  • \( \Pr(\text{Mam=+} | \text{Cancer=+}) = 0.90 \)
  • \( \Pr(\text{Mam=+} | \text{Cancer=-}) = 0.07 \)

• Maximum de vraisemblance : quelle est la valeur de Cancer qui maximise \( Mam=+ \) ?

  • \( \Pr(\text{Mam=+} | \text{Cancer=+}) = 0.90 \)
  • \( \Pr(\text{Mam=+} | \text{Cancer=-}) = 0.07 \)

• Considérons 1000 femmes âgées de 40 à 50 ans, sans antécédents familiaux et sans symptômes de cancer
  • 8 femmes sur 1000 ont un cancer
• On réalise une mammographie
  • Sur les 8 femmes ayant un cancer, 7 auront un résultat positif
  • Sur les 992 femmes restantes, 69 auront un résultat positif
• Une femme passe une mammographie, le résultat est positif
• Que devrait-on inférer ?

\[ \Pr(\text{Cancer=+} | \text{Mam=+}) = \frac{7}{7 + 69} = \frac{7}{76} \approx 0.09 \]

\[ \color{purple}{p(\theta \vert x)} = \dfrac{\color{orangered}{p(x\vert \theta)} \color{steelblue}{p(\theta)}}{\color{green}{p(x)}} \]

\( \color{steelblue}{\Pr(\theta)} \) la probabilité a priori de \( \theta \) : tout ce qu'on sait de \( \theta \) avant d'observer les données. En l'occurrence : \( \Pr(\text{Cancer=+}) = 0.008 \) et \( \Pr(\text{Cancer=-}) = 0.992 \).

prior <- c(0.008, 0.992)

\[ \color{purple}{p(\theta \vert x)} = \dfrac{\color{orangered}{p(x\vert \theta)} \color{steelblue}{p(\theta)}}{\color{green}{p(x)}} \]

\( \color{orangered}{p(x\vert \theta)} \) probabilité conditionnelle des données \( x \) sachant le paramètre \( \theta \), qu'on appelle aussi la likelihood (ou fonction de vraisemblance) du paramètre \( \theta \).

like <- rbind(c(0.9, 0.1), c(0.07, 0.93) ) %>% data.frame
colnames(like) <- c("Mam+", "Mam-")
rownames(like) <- c("Cancer+", "Cancer-")
like

        Mam+ Mam-
Cancer+ 0.90 0.10
Cancer- 0.07 0.93

\[ \color{purple}{p(\theta \vert x)} = \dfrac{\color{orangered}{p(x\vert \theta)} \color{steelblue}{p(\theta)}}{\color{green}{p(x)}} \]

\( p(x) \) la probabilité marginale de \( x \) (sur \( \theta \)). Sert à normaliser la distribution.

\[ \color{green}{p(x)=\sum\limits_{\theta}p(x|\theta)p(\theta)} \]

(marginal <- sum(like$"Mam+" * prior) )

[1] 0.07664

\[ \color{purple}{p(\theta \vert x)} = \dfrac{\color{orangered}{p(x\vert \theta)} \color{steelblue}{p(\theta)}}{\color{green}{p(x)}} \]

\( \color{purple}{p(\theta \vert x)} \) la probabilité a posteriori de \( \theta \) sachant \( x \), c'est à dire ce qu'on sait de \( \theta \) après avoir pris connaissance de \( x \).

(posterior <- (like$"Mam+" * prior ) / marginal )

[1] 0.09394572 0.90605428

Avant de passer le mammogramme, la probabilité qu'une femme tirée au sort ait un cancer du sein était de \( \Pr(\text{Cancer=+}) = 0.008 \) (prior). Après un résultat positif, cette probabilité est devenue \( \Pr(\text{Cancer=+} | \text{Mam=+}) = 0.09 \) (posterior). Ces probabilités sont des expressions de nos connaissances. Après un mammogramme positif, on pense toujours que c'est “très improbable” d'avoir un cancer, mais cette probabilité a considérablement évolué relativement à “avant le test”.

A Bayesianly justifiable analysis is one that treats known values as observed values of random variables, treats unknown values as unobserved random variables, and calculates the conditional distribution of unknowns given knowns and model specifications using Bayes’ theorem (Rubin, 1984).

Il s'agit d'un problème de probabilités conditionnelles… Définissons les événements suivants :

P1 : l'animateur ouvre la porte 1
P2 : l'animateur ouvre la porte 2
P3 : l'animateur ouvre la porte 3

V1 : la voiture se trouve derrière la porte 1
V2 : la voiture se trouve derrière la porte 2
V3 : la voiture se trouve derrière la porte 3

Si on a choisi la porte n°1 et que l'animateur a choisi la porte n°3 (et qu'il sait où se trouve la voiture), il s'ensuit que :

\( \Pr(\text{P3} | \text{V1}) = \dfrac{1}{2}, \quad \Pr(\text{P3} | \text{V2}) = 1, \quad \Pr(\text{P3} | \text{V3}) = 0. \)

On sait que \( \Pr(\text{V3} | \text{P3}) = 0 \), on veut connaître \( \Pr(\text{V1} | \text{P3}) \) et \( \Pr(\text{V2} | \text{P3}) \) afin de pouvoir choisir. Résolution par le théorème de Bayes.

\( \Pr(\text{V1} | \text{P3}) = \dfrac{\Pr(\text{P3} | \text{V1}) \times \Pr(\text{V1})}{\Pr(\text{P3})} = \dfrac{\dfrac{1}{2} \times \dfrac{1}{3}}{\dfrac{1}{2}} = \dfrac{1}{3} \)

\( \Pr(\text{V2} | \text{P3}) = \dfrac{\Pr(\text{P3} | \text{V2}) \times \Pr(\text{V2})}{\Pr(\text{P3})} = \dfrac{1 \times \dfrac{1}{3}}{\dfrac{1}{2}} = \dfrac{2}{3} \)

Nos intuitions probabilistes sont, dans la grande majorité des cas, très mauvaises. Au lieu de compter sur elles, il est plus sage de se reposer sur des règles logiques (modus ponens et modus tollens) et probabilistes simples (règle du produit, règle de la somme, théorème de Bayes), nous assurant de réaliser l'inférence logique la plus juste.

Autrement dit, “Don't be clever” (McElreath, 2015).